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Abstract 

The energy distributions of secondary electrons produced by energetic carbon ions (in the energy 
range used, e.g., in hadron therapy) , incident on liquid water, are discussed. For low-energy ions, 
a new parameterization of the singly-differential ionization cross sections is introduced, based on 
tuning the position of the Bragg peak. The resulting parameterization allows a fast calculation of 
the energy spectra of secondary electrons at different depths along the ion's trajectory, especially 
near the Bragg peak. At the same time, this parameterization provides penetration depths for 
a broad range of initial-ion energies within the therapeutically-accepted error. For high-energy 
ions, the energy distribution is obtained with a use of the dielectric response function approach. 
Different models are compared and discussed. 
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I. INTRODUCTION 



An accurate energy distribution of electrons produced as a result of the ionization of water 
molecules induced by fast incident ions is a key piece of information for many applications. 
One of the most important among these applications, since liquid water is considered to 
be a good tissue-like medium, is related to its relevance to the physics of ion-beam cancer 
therapy and radiation protection in space. Ion beam cancer therapy has become an operative 
treatment tool, whose power is based on theoretical and computational methods of basic 
science On the other hand, radiation damage by ions is becoming a major topic of 
research for shielding of human space missions by the ESA and NASA [2j. These two fields 
of applications differ slightly in the initial energies of incident ions, hundreds of MeV/u for 
cancer therapy and ~1 GeV/u for radiation protection from galactic cosmic rays. 

The goal of studies applied to both curing and shielding purposes is an understanding and 
controlling of the radiation damage, which is focused on DNA damage and its repair. Many 
details of the involved physical processes, starting from the incidence of an energetic ion on 
tissue and leading to biological damage, are still far from being explained and calculated on 
a nanoscopic level. An attempt to build an inclusive picture of relevant physical processes 
was suggested in our previous work [3|, |4j, which proposed a multiscale approach to this 
problem. 

It is widely accepted that the major part of damage done by ions is related, directly or 
indirectly, to the secondary electrons (frequently called (^-electrons) produced by ionization 
of the medium. These electrons may interact with parts of DNA molecules in the cell nuclei, 
generate other secondaries, such as other electrons or radicals, which can then interact with 
DNA. Besides this, electrons carrying energy, with their propagation define a volume around 
ihe track of abrupt temperature increase that can lead to a local heating induced damage 
5)]; furthermore each ionization produces a hole which may also cause DNA damage. 

In the multiscale approach {3), the energy spectrum of secondary electrons occupies the 
key position. It is also an important input for numerous Monte Carlo (MC) simulations of 



track structure 



7 



] . An increasing interest for a more accurate shape of this distribution, 
especially at low energies, has followed the well known experiments by the Sanche's group [s| 
revealing possible damage induced by low-energy (~5 eV) electrons. 

The studies of the energy distribution of secondary electrons have been carried out by 
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different groups (e.g. [a. llOip . In a majority of works, energetic electrons or protons are the 



primary projectiles |10l. Illl]; while a few studies have been extended to heavier ions 12|, [13] 
In our analysis, we consider carbon ions because they have been largely used in ion-beam 
therapy both in Germany and Japan [141 ] . However, our analysis can also be extended to 
other ions. 

The difficulties in analyzing singly differential ionization cross sections are basically of 
two types: (1) the impossibility to depict all energy ranges using the same formalism be- 
cause the energy range of interest is too large to be included in the range of validity of a 
single approximation; (2) the deviation of the ion's propagation through liquid water from 
that of water vapor often used in experiments, is unknown. As the applicability of the Born 
approximation is a key issue in this connection, we divide our analysis (as it is usual in 
similar works [si, [h]]) between the fast and relatively slow regimes of the projectile's veloc- 
ity, corresponding to different parts of the ion's trajectory. First, we start with a general 
approach valid for slow ions. This approach is the only possibility in that energy range, and 
the accuracy in the shape of secondary electron spectra is limited in the entire energy range. 
Then, we present a method valid only for fast ions, but returning more accurate electron 
energy distributions in that range. 

II. SLOW REGIME 

In order to account for the region where the ion's energy is smaller than 500 keV/u, 
where the Born approximation is not valid, we use a parametric approach built on existing 



experimental data. As was shown in our previous papers 15j, ]l6| , the probability to produce 



dN secondary electrons with kinetic energies from W to W + dW, j^dW, from a segment 
of the ion track d(, is proportional to ^d(, where ^ is the singly differential ionization 
cross section (SDCS). This cross section is the main characteristic in our analysis defining 
both the ion's stopping in the medium and the energy spectrum of secondary electrons. 
Since the angular distribution of prevailing low-energy electrons (W < 45 eV) is rather 



flat over the whole angular range 171], the SDCS is sufficient for our analysis. Besides 
the kinetic energy of secondary electrons and the properties of water molecules, the SDCS 
depends on the kinetic energy of projectiles, T, and their charges, z. 
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TABLE I: Fitting parameters for the two inner shells of the water molecule (Icq, 2ai ) an d for the 
three outer shells (lb 2 , Sa±, lb\) in the original parameterization by Rudd for vapor 
new parameterization here introduced for liquid. 



18l and in a 



parameter 


A 1 B 1 d D 1 E l A 2 B 2 C 2 D 2 a 


inner shell 


1.25 0.5 1.0 1.0 3.0 1.1 1.3 1.0 0.0 0.66 


outer shell vapor [18] 


0.97 82 0.4 -0.3 0.38 1.04 17.3 0.76 0.04 0.64 


outer shell liquid (this work) 


1.02 82 0.50 -0.78 0.38 1.07 14.5 0.61 0.04 0.64 



A. Semi-empirical Model for vapor water 



In our previous papers 15j, [l6| , we used a model based on a parameterization by Rudd 18] , 
which yields an analytical formula for SDCS on a wide energy range combining the exper- 
imental data for water vapor, calculations within the plane-wave Born approximation, and 
other theoretical models. This model, original for protons, has been extended first to he- 



lium [12| and then to heavier ions by us 
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We used the following parametric functions 

da(W } T) 



and a few other works 



16j: 



19] 



dW "eff ^ j 

1 L 



(1) 



(1 + Uif (1 + exp(c^ - < ax )M)) 



F{(v) = A\ 



Hi 



l+v 1 



P 2 



+ 



B\/v 2 + v 2 1 + E[v D \ +A 



F l 2 {v) 



A\v 2 + B\ 



C> D 2+4 + Aiv 2 + Bi 



(2) 
(3) 



where the sum is taken over the electron shells of the water molecule, a$ is the Bohr radius, 
R is the Rydberg, iV, the shell occupancy, Ii the ionization potential of the shell, Ui = Wjl% , 
Vi = (mV 2 / (2/j)) 1 / 2 , cj 4 max = 4v 2 — 2vi — Rj (4/j); m is the mass of electron, V is the velocity 
of the projectile, T its kinetic energy and (3 = V/c. 

The corresponding fitting parameters, taken from Ref. [3], A\ ... E{, A\ ... D 2 , a 1 are 
listed in Table [Hand the ionization potentials in Table ILT1 fv apor row). 

This parameterization is different from the original [18( in two respects: first, the ex- 
pression for F\ is modified so that it has the correct asymptotic behavior in the relativistic 
limit [3], and second, it contains an effective charge of the ion z e ff, which depends on the 
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TABLE II: Ionization potentials for vapor 18| and liquid [3j] water (in eV). In the two approaches 
presented in this paper, the first and the second set are respectively used. 



shell 


Ibi 


3cq 


1&2 


2ai 


lai 


vapor 


12.61 


14.73 


18.55 


32.2 


539.7 


liquid 


10.79 


13.39 


16.05 


32.3 


539.0 



velocity of the ion. The effective charge takes into account the effect of charge transfer, i.e., 
gradual reduction of the original charge of the ion as it slows down. For the effective charge, 

z(l -exp(-125/3z- 2/3 )). 



we have used an expression given by Barkas 



20] 



Z eff 



These two modifications give the shape of the LET dependence along the track for a single 
ion. However, ions in the beam are spread in energies due to scattering. The straggling in 



energies is taken into account in Ref. [l6j, but more details are given in Ref. [2l| . All in 
all, our resulting general parameterization of SDCS for water yielded a 97% accuracy in 
reproducing the Bragg's peak position in comparison with experiment and compared very 
well with MC simulations 22j in shape (without nuclear fragmentation). 



B. Semi-empirical model extensions for liquid water 



However, the parameterization for SDCS, used in Refs. [16j, [2l|], did not take into ac- 
count medium density effects, i.e., the effects arising from the difference between the ion's 
propagation in liquid water from that of vapor. We used the liquid-water density, but took 
all other parameters corresponding to water vapor. Unfortunately, the experimental data 
for ionization cross sections in liquid water are very limited and insufficient for making an 
independent parameterization. However, since the SDCS of ionization along with excita- 
tion contribute to the stopping cross section, directly related to LET 16|, it is possible to 



extract a modified parameterization of the SDCS for liquid water from different available 



experiments. In Ref. [k]], the experimental measurements of stopping cross sections 23] 
are used as a function of the ion's energy. This only provides a single curve to which to 
compare. Instead, we choose to fit the parameterization with measurements of the LET 
dependence on the penetration depth in liquid water. These depths have been measured 



at GSI 



24 



251 ] . This not only provides several curves at different initial energies, but also 



a direct comparison with the carbon ion data. So instead of tuning parameters directly to 
the stopping cross sections as in n|, we use the reported Bragg peak profiles, for different 
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TABLE III: Energies and parameters for excitation cross sections (Eq.[6|) corresponding respectively 
to states A x Bi, B 1 Ai, Rydberg A+B, Rydberg C+D, and diffuse bands. 



k 


E k (eV) 


O-k bk Ck 


1 


8.17 


2245 4493 0.85 


2 


10.13 


6319 7020 0.88 


3 


11.31 


4387 8135 0.88 


4 


12.91 


989 3172 0.78 


5 


14.50 


1214 3352 0.78 



initial ion energies and thus the LET as a function of the penetration depth. 

In fact, the integration of the inverse of the LET gives the position of the Bragg peak. 
Hence, we can adjust the parameters of the SDCS by tuning the depth of the Bragg peak 
to the experiment. Such a procedure is rather simple and fast using our approach 16]. 

So, we start with the SDCS of ionization, which includes relativistic effects and charge 
transfer as an input to the stopping cross section, defined as 



a. 



St 



E 



(4) 



Now, we need to add the contribution of excitation to the stopping cross section. These cross 
sections are accounted by scaling with our effective charge, depending on ions' velocities, a 
semi-empirical formulation by Miller and Green 26| for protons 



@exc,k 



2 a k (T — E k 

Zeff b k + T c * 



(5) 



where the index k corresponds to five different excitation transitions, i.e., to states A 1 B 1 , 
B 1 A i , Rydberg A+B, Rydberg C+D, and diffuse bands, respectively; a exCyk {T) is the excita- 
tion cross section, E exCyk is the corresponding excitation energy and a k ,b k , c k are parameters 
reported in Table III and computed from data reported and fully explained in Refs. [lOj, |26j , 



giving cross sections in A 2 when energy is input in eV. The stopping cross section is then 
calculated [16| as 



(T* t (T) = CT st (T) + ^ E exc,i&exc,i( T ) ■ ( 6 ) 
i 

The LET is given by na* t (T), where n is a number density of liquid water molecules. 
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FIG. 1: (Color online) Linear energy deposition (left panel) for a carbon ion with incident energy 
of 400 MeV/u, computed with different parameterizations, within the experimental error of the 
original data, and compared to experiments [25| and MC simulations [zj (where no nuclear frag- 
mentation is considered). The zoom in the inset shows also the experimental penetration depth. 
On the right panel the effect of this new parameterization (solid line) on the spectra of secondary 
electrons compared to the previous one (dashed line). 



At this point, we add the energy straggling correction, i.e., a widening of the peak arising 
from the stocasticity of the energy loss phenomenon and related to the multiple scattering of 
ions. This correction affects also the position of the Bragg peak, slightly reducing its depth. 
The correction is given by 



dx 



(7) 



cr str (x )v / 27T 



■I'D 



f (x')ex P ( 



(x 1 - x) 2 



2a str (x y 



)dx' 



where xq is a maximum penetration depth of the projectile, and <J s tr(xo) — 0.012xg 951 /-\/]4 
is the longitudinal-straggling standard deviation (in cm, when xq is also in cm) computed 
with a phenomenological formula introduced by Chu [27|, and depending (through x ) on 
the initial energy of the ion To. A is the ion's mass number. We found, that if we gradually 
change the parameters from Rudd's model within the experimental error ljj] , the Bragg peak 
position can be tuned to match those obtained in stopping power experiments (see Fi g, [p , 
as well as MC transport code simulation (MCHIT code based on GEANT4 toolkit 0, 28]). 
The tuned parameters are shown in the third row of Table [B The same parameterization 
is applicable to a broad range of energies (see Table IIVI and Fig. [2]), keeping the Bragg 
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FIG. 2: Linear energy deposition for carbon ions for different initial energies To, with our model 
compared to experiments from GSI 24|, [25]. 



TABLE IV: Bragg peak positions for different incident energies of carbon ions in water obtained 
with the present parameterization with and without accounting for straggling and compared with 
experiments 2J, |25|] and Monte Carlo simulations 22 



2$. 



T(MeV/u) 


Bragg peak position (mm) 




us (single ion) 


us (straggling) 


Schardt 


MCHIT 


400 


275.75 


275.01 


274.72 


274.53 


330 


201.81 


201.30 


201.42 


200.52 


270 


144.78 


144.42 


144.82 


143.78 


200 


87.19 


86.94 


86.46 


86.90 


195 


83.49 


83.26 


83.39 




135 


44.20 


44.05 


43.34 


43.73 



peak position within the acceptable error for therapeutic uses (0.5 mm). In particular, the 
comparison with Monte Carlo simulations, shows on one hand an almost perfect agreement 
in the shape of the curves when considering only the electromagnetic processes and no 
nuclear fragmentation channels; as shown in the limit case of a 400 MeV/u ion beam (left 
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panel of Fig. [TJ), where fragmentation has a drastic effect on the peak shape, but nevertheless 
doesn't affect the depth. On the other hand it reveals a slightly better reproducibility of the 
experimental depths in the whole range of energies for the present method (Table HV1) . as 
compared to MC, which is based on a single value of average excitation energy (< / >=85 
eV). 

Thus, the approach suggested in Ref. [16J with the new parameters given in Table [H 
can be used successfully also for rapid calculations of the Bragg peak position; and, in 
case of limited nuclear fragmentations, as occur in smaller ions, may also provide a correct 
deposition profile depiction [29] 



III. FAST REGIME 



When ions are more energetic than 500 keV/u, the Born approximation is valid, and 
a direct connection with the experimental data for liquid water is possible through the 
dielectric response model. This provides an opportunity for a more accurate calculation of 
the SDCS and thus the energy spectrum of secondary electrons. 



A. The dielectric response model: optical approximation 



In this model, presented by a number of studies for electrons and protons [10j, |3Cj , but 
rarely extended to heavier ions 13J , the only parametric quantity is retrieved from measure- 
ments of photoionization cross sections. The main advantage of this model, which makes it 
suitable, especially for treating condensed media, is the simultaneous accounting for both 
single-particle and collective effects in the analysis of the response of a medium to the ioniz- 
ing/exciting energetic particle. The key quantity is the energy loss function, often denoted 



as 7/ 2 , derived by the die 



strength df/dE(E,k) 
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ectric function e and directly connected to the generalized oscillator 



22| : 



rj 2 (E, k) = Im 



-1 



e{E,k) 



nn 2 p df(E,k) 
2ZE dE 



(8) 



where Z is the total number of electrons of a molecule of the medium, E is the transfered 



energy, k is the transfered momentum, and VL p = ARina^nZ) 1 / 2 = 21.46 eV is the plasma 
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frequency computed for liquid water where Z = 10. 

Then, for a liquid medium, it is more convenient to insert a macroscopic cross section 



A = ncr, and consider 



32] 



dA _ z 2 eff m p7]2 (E,k) 



(9) 



dEdk 7ra T k 

where the factor m p is the mass ratio between proton and electron, and, after integration 
over k 



dA Z eff m p 



kmax 



, ^ i \ d/kj 
r]2{E,k) — 



dE ITCLnT lu 

This can be expanded under the condition E « T/m p in the following sum 



(10) 



dA 



2 

Z eff m P 



dE ttclqT 



V -^^ln^- + B{E) + 0(Em p /T) 



(11) 



In general, e(i£, fc) = e±(E, k) + iei{E, k), where both the real part e\ and the imaginary 



part €2 are derivable from experiments 



33( . They are shown in Fig. [3] in the optical limit 



(k=0). Based on these functions, it is possible to derive r]2 vital for our calculations. 

Thus, in cases where the optical limit is reasonable, i.e., when it is possible to consider 
ionization by ions to be similar to that of photons, it is possible to get a good approximation 
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of the macroscopic cross section by retaining only the first term of Eq. (fTTjl . which is called 
Bethe optical term 

In passing from the energy loss to the secondary electron distribution, it is necessary to 
introduce the different ionization thresholds Ij (here the real liquid water data are used, 
second row of Table [III) an d contributions from different shells 

dA ^ dAj 
2-^ 7w 



dW ^ dW 

j 

dA j „ z hf m p 



dW na T 



2 m p R 



(12) 



The factors Gj (E) account for different contributions to the optical spectrum from ionization 
of different shells. They are computed as suggested in a series of studies performed for 
protons Ref. 30| by treating separately the different transitions by decomposing the resulting 



optical data in a sum of Drude functions, with the constrainth of respecting the sum rule. 

ion exc 

e 2 (E, 0) = 4!L(£, 0) + Y, €*c(E, 0) (13) 

j i 

where the e^] on and the e^ exc represent, respectively, the contributions of a given ionization 
or excitation mode to the total imaginary part of the dielectric function and are computed 
differently. In particular the ionization contributions are expressed well by conventional 
Drude functions: 



where the parameters fj,jj, Ej are respectively oscillator strength, width and position of 
the given transition and are found by imposing the sum rules: 

POO 

/ Ee 2 (E,0)dE = n 2 p n/2 (15) 
Jo 



POO 

/ E V2 (E,0)dE = n 2 p n/2 (16) 
Jo 
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TABLE V: Parameters ( pos itions, widths and strengths) for the ionization transitions used for the 
shell contributions to €2 



shell 


^■(eV) 7i (eV) 


fi 


I61 


16.30 


14.00 


0.2300 


3a± 


17.25 


10.91 


0.1600 


I62 


28.00 


27.38 


0.1890 


2oi 


42.00 


28.68 


0.2095 


lai 


450 


360 


0.3143 
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FIG. 4: Energy spectra of secondary electrons induced by a carbon ion: Comparison between the 
parametric model and different dielectric response model approaches for liquid water. 



We use here the last fitted parameters from Ref. [llj (see table 5), for the single contri- 
butions, while we directly derived ex and e 2 from experiments in the region 1-50 eV, and 
extrapolating with a non-linear curve fitting to the asymptotes (respectively 1 and 0). Thus 
we take Gj(E) = e^l^E , 0) / e 2 (E , 0) . We can probe this approximation, as mentioned, un- 
der the condition of E <C T/m p ; this turns out to be valid for relatively high T (>0.5 MeV), 
when considering the lower energy part of the secondary electron spectrum. From Fig. HI 
we can see that compared to our parametric formalism, there is a region where the two ap- 
proaches give very close results (10-40 eV), while at a lower energy the curve corresponding 
to the parameterization deviates, as it is not able to correctly reproduce the liquid medium, 
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FIG. 5: Energy spectra of secondary electrons induced by a carbon ion: Comparison between the 
reported models and Ref. fl3j ]. 



and at a higher electron energy the Bethe optical approximation is not valid; thus, the 
correct behavior is there described by the parametric, vapor-like curve. 



B. The dielectric response model: extension in the energy-momentum plane 

Finer approaches extend to k >0, i.e., gen erating the full energy-momentum plane, char- 



acterized by the well-known Bethe ridge 



341 ]. through different models, usually called dis- 



persion schemes. These models follow the requirements to generate a generalized energy 
loss function 1] 2 (E, k) obviously consistent in the optical limit (k — > 0) to the corresponding 
optical function, as well as in the asymptotic behaviour for large momentum transfer to the 
free electron limit: 

lim r) 2 (E, k) ~ 5(E - Q) , (17) 

k— >oo 



where Q = k 2 /2m is the free electron recoil energy 35]. Many methods have been proposed, 
such as the one by Ashley [3] , based on a 5-function representation of the optical oscillator 
strength of an atom 
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ViE> = Y,fAE-E i ) (18) 

j 

(with fj and Ej respectively height and position of a given energy transition), where Ej is 
extended for k > to Ej + Q and a generalization to condensed matter is performed: 

/OO 771/ 771 ,0 

-772^', 0)5(£ - (£' + Q))dE' = — - Q, 0)Q(E - Q) (19) 

where Q(E — Q) is a step function, equal to 1 for E > Q and otherwise. The k dependence 
from the optical data is then achieved. 

Assuming that the weights computed in the case of momentum transfer don't change 
significantly, it is then possible to assume G(E, k) = G(E) and compute the contributions 
to ionization of the shells with the eq.19 and 

dAj(E,k) z lff m p f kmax „ ,„> ,-, ..dk 



?//' 

dW 7ia T 



/ Gj(E)r] 2 (E, k)— . (20) 

v krn.i.n. 



The integration in k, as shown in Eq. fflQl) . is then performed with the minimum and 
maximum momentum transfer as the integration limits: 

k± = v / 2^(v / T ± y/T-E) (21) 

Thorough reviews of the various possibility to expand in the Bethe surface are available 



in Refs. 



37 



Other approaches include also a way for accounting the relativistic regime 
of the projectile energy. In general, this correction is small in the energy ranges where 
is the main focus of interest in evaluation of electrons spectra (Bragg peak region, i.e. 
T < lOMeV/u), where f3 2 < 10~ 2 . We report here, anyway, the results obtained by using 



such an approach, suggested by Asoskov et al. |39j], extended by us from an atomic single 
shell system to include a multishell target. This formalism uses a different method but 
perform still a similar expansion in the energy-momentum plane by accounting for the recoil 
energy in order to express the generalized oscillator strength (GOS) as a function of the 
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optical oscillator strength (OOS). That, after integration in k, bring to the formula: 



dA 
dW 



dAj 
dW 



EdAj 
dW 



8a 2 R 2 N t z eff Z (Ee^(E) 



VL 2 p mV 2 



Ee 2 



In 



2mV 2 



2^2 



E\l ~(3 2 e\ 



£2 



■ arctan 



-P 2 e 2 



1 - /? 2 ei 



(22) 



+ 



F{E ) 
E 2 



where 



F(E) 



E E'e^UE' 



■dE' 



(23) 



and, for each shell, E = W + Ii- 

In Fig. HI one can appreciate corrections induced by these approaches to the Bethe optical 
term, which is reliable only for very small values of W, as expected. We can see how the 
parametric and other curves are almost coincident at large W, where the medium differences 
are less dramatic. Finally, in Fig. [5] we report a comparison with another model 13]. 



Ref. 



13) uses a stochastic simulation for the partition of the different shells, and report only 



normalized values (hereby scaled for the present total cross section). In this case the results 
are slightly The major difference between the results is the slightly lower values obtained 
in Ref. 13j at very low energy. This region of course is still very important, thus a further 
improvement of detail of the energy shape, with the help of further experiments, is desirable. 



IV. CONCLUSIONS 



This paper gives an opportunity to calculate the energy spectra of secondary electrons 
produced by energetic ions incident on tissue, mimiced by liquid water. In order to report 
these electron energy distributions with the best possible accuracy for each ion energy range, 
we consider different approaches. 

The energy distributions of secondary electrons produced by low-energy ions are obtained 
using a general parameterization based on the Rudd Model with parameters modified for 
liquid water distinguished from the original parameters fitted to experiments with water 
vapor. The new parameters were obtained by fitting the depths of Bragg peaks (for different 
incident energies) to the experiments done at GSI. The obtained parameterization can be 
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used for a fast calculation of the position of the Bragg peak for a given energy of the 
projectile with a therapeutically accepted precision. An advantage of such parameterization 
is its universality for different applications and its analyticity, which makes all calculations 
fast. 

For higher ion energies, the dielectric response approach, previously tested for protons, 
and herein applied for heavier ions, is shown to be successful in describing the details of 
secondary electron spectra profiles. We compare two approaches and show that they can be 
used for different values of secondary electron energies. 

The results of this work can be used for calculations of realistic radial dose distributions at 
different depths along the ion's trajectory. Such calculations, if rapid enough, are important 
for calculating the relative biological effectiveness and treatment planning. 

Further experiments, combining more realistic projectiles (heavier ions, since the direct 
ionization data are only available for protons and helium ions 



171 ]) and more realistic media, 



are very much desired. 
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